##__docformat__ = "restructuredtext en"
# ******NOTICE***************
# two_opt.py module by Shoibal Chakravarty
#
# You may copy and use this module as you see fit with no
# guarantee implied provided you keep this notice in all copies.
# *****END NOTICE************
"""
This module  contains the functions that calculate the utility functions for the scenarios
with regions that have one period look ahead planning horizon.
"""

__author__ = "Shoibal Chakravarty (shoibalc@princeton.edu)"
__version__ = "0.1"
__date__ = "2011/August"
__license__ = "BSD"
__all__ = ['U2','U2optimized','U3','U3optimized']
import numpy as np
from dataops import cg2


def U2(gx0,*args):
    """
    U2 calculates the two period utility of region x for a given value of
    region y's green deployment, region c has a myopic strategy that is known. U2 is used to calculate 
    the best response of x (i.e.  deployment with max utility) for a given y.\n\n
    Function parameters:\n\n
    t = current period (t+1 = next period)\n
    x = region with a Nash strategy whose best response is to be determined\n
    y = the other region with a Nash strategy\n
    c = third region with a myopic strategy\n\n
    gx0 = new green energy of region x in period t \n
    gy0 = new green energy of region y in period t \n
    Bx = green preference of region x in period t\n
    Bx1 = green preference of region x in period t+1\n
    cnx = levelized cost of new fossil in region x in period t\n
    cnx1 = levelized cost of new fossil in region x in period t+1\n
    cox = levelized cost of old fossil in region x in period t\n
    cox1 = levelized cost of old fossil in region x in period t\n
    Gx0 = total green energy in region x in period t\n
    Gy0 = total green energy in region y in period t\n
    Gc0 = total green energy in region c in period t\n
    Gc1 = total green energy in region c in period t+1\n
    G0 = total global green energy initially\n
    ex0 = new energy in region x in period t\n
    Ex0 = total energy in region x in period t\n
    Ex1 = total energy in region x in period t+1\n
    M: min( M*Ex0,Ex0 - Gx0) is an upper bound on gx0
    M1: min( M1*Ex1,Ex1 - Gx0 - gx0) is an upper bound on new green energy in region x in period t+1\n
    delta = discount factor of period t+1 wrt period t\n
    detail = boolean flag controlling function return parameters
    """
    gy0,cg0,gpar,Bx,Bx1,cnx,cox,cnx1,cox1,Gx0,Gy0,Gc0,Gc1,G0,ex0,Ex0,ex,Ex1,M,M1,delta,detail = args

#    mt0 = max(min(M*Ex0,(Ex0 - Gx0)),ex0)
    mt0 = min(M*Ex0,(Ex0 - Gx0))
    ux0, cx0 =  0.0, 0.0
    eps = 1e-8
    Gt = Gx0 + Gy0 + Gc0
    cg = cg2(cg0,G0,Gt,gpar)
    cgn = cg - cnx
    cgo = cg - cox
    if cgn > 0.0:
        if mt0 < ex0-eps:
            print mt0,ex0,Gt,Ex0,Gy0,Gx0
            raise ValueError('Max deployment should be greater than next period net demand. Check your assumptions!\n')
        if gx0 <= ex0:
            cx0 = gx0*cgn
            ux0 = Bx*np.log(1.0 + gx0) - cx0 
        elif gx0 > ex0:
            cx0 = (gx0-ex0)*cgo + ex0*cgn
            ux0 = Bx*np.log(1.0 + gx0) - cx0
        else:
            print  gx0, ex0, mt0
            raise ValueError("function U2 1st period: part 1 cgn > 0.0, something wrong! you should never reach here!\n")
    elif cgn <= 0.0 and cgo > 0.0:
        if gx0 <= ex0:
            cx0 = 0.0
            ux0 = Bx*np.log(1.0 + gx0)
        elif gx0 > ex0:
            cx0 = (gx0-ex0)*cgo
            ux0 = Bx*np.log(1.0 + gx0) - cx0
        else:
            print  gx0, ex0, mt0
            raise ValueError("function U2 1st period: part 1 cgo > 0.0, something wrong! you should never reach here!\n")
    elif cgo <= 0.0:
#### This hack prevents the "consumption smoothening"" type behavior that we see in some cases. After we reach co < cg we deploy
#   the maximum possible in that period, i.e. the model effectively becomes "myopic".
        if gx0 >= mt0 - 0.5:
            cx0 = 0.0
            ux0 = Bx*np.log(1.0 + gx0)
#            print "SPECIAL CASE", gx0,mt0,ux0
        else:
            cx0 = 0.0
            ux0 = -1e8
#        cx0 = 0.0
#        ux0 = Bx*np.log(1.0 + gx0)
#############################################
    else:
#        print  gx0, cg,cox0, cnx0, ex0, mt0
        raise ValueError("function U2 1st period: part 2, something wrong! you should never reach here!\n")
    
#    print "g", gx0, "u", ux0, "e", ex0, "m", mt
    Gx1 = Gx0 + gx0
    Gy1 = Gy0 + gy0
#    mt = max(min(M1*Ex1,(Ex1 - Gx1)),ex)
    mt = min(M1*Ex1,(Ex1 - Gx1))
    Gt = Gx1 + Gy1 + Gc1
###    cg = cgx*(Gt/G0)**alpha
    cg = cg2(cg0,G0,Gt,gpar)
    cgn = cg - cnx1
    cgo = cg - cox1
# Change B to next period value
    Bx = Bx1
    g1  = Bx/cgn - 1.0
    g2  = Bx/cgo - 1.0
    
    if mt < ex-0.00001:
#        print mt,ex,Gt,Ex,Gy,Gx
        raise ValueError('Max deployment should be greater than next period net demand. Check your assumptions!\n')
    if Bx >= cgn:
        if cgn > 0.0:
            if g1 < ex:
                greturn = g1
                creturn = greturn*cgn 
                ureturn = Bx*np.log(1.0 + greturn) - creturn
            elif g1 > ex and g2 < ex:
                greturn = ex 
                creturn = greturn*cgn
                ureturn = Bx*np.log(1.0 + greturn) - creturn
            elif g2 > ex and g2 <= mt:
                greturn = g2
                creturn = (greturn-ex)*cgo + ex*cgn
                ureturn = Bx*np.log(1.0 + greturn) - creturn
            elif g2 > mt:
                greturn = mt
                creturn = (greturn-ex)*cgo + ex*cgn
                ureturn = Bx*np.log(1.0 + greturn) - creturn
            else:
                raise ValueError("U2 part 2 something wrong! cgn > 0.0: you should never reach here!\n")
#                raise ValueError("nextgG something wrong! cgn > 0.0: you should never reach here!\n")
#            print "cg", cg, "g", greturn, "u", ureturn, "c", creturn, "g1",g1, "g2", g2,"e", ex, "m", mt
#            return greturn, ureturn, creturn, mt
        elif cgn <= 0.0 and cgo > 0.0:
            if  g2 < ex:
                greturn = ex
                creturn = 0.0
                ureturn = Bx*np.log(1.0 + greturn)
            elif g2 > ex and g2 <= mt:
                greturn = g2
                creturn = (greturn-ex)*cgo
                ureturn = Bx*np.log(1.0 + greturn) - creturn
            elif g2 > mt:
                greturn = mt
                creturn = (greturn-ex)*cgo
                ureturn = Bx*np.log(1.0 + greturn) - creturn
            else:
                raise ValueError("U2 part 2 something wrong! cgo > 0.0: you should never reach here!\n")
#                raise ValueError("nextgG something wrong! cgo > 0.0: you should never reach here!\n")
#            print "cg", cg, "g", greturn, "u", ureturn, "c", creturn, "g1",g1, "g2", g2,"e", ex, "m", mt
#            return greturn, ureturn, creturn, mt
        elif cgo <= 0.0:
            greturn = mt
            creturn = 0.0
            ureturn = Bx*np.log(1.0 + greturn)
#            print "cg", cg, "g", greturn, "u", ureturn, "c", creturn, "g1",g1, "g2", g2,"e", ex, "m", mt
#            return greturn, ureturn, creturn, mt
        else:
            raise ValueError("U2 part 2 something wrong! Bx > cgn: you should never reach here!\n")
    elif Bx < cgn:
        greturn = 0.0
        ureturn = 0.0
        creturn = 0.0
#        print "cg", cg, "g", greturn, "u", ureturn, "c", creturn, "g1", g1, "g2", g2,"e", ex, "m", mt
#        return greturn, ureturn, creturn, mt
    else:
        raise ValueError("U2 part 2 something wrong! you should never reach here!\n")
    if detail:
#        print "g", greturn, "u", ureturn, "g1",g1, "g2", g2,"e", ex1, "m", mt
        return  delta*ureturn + ux0, gx0, ux0, cx0, greturn, ureturn, creturn, mt0
    else :
#        print "g", greturn, "u", ureturn, "g1",g1, "g2", g2,"e", ex1, "m", mt,"g0", gx0, "u0", ux0
#        print "gx0", gx0, "g2", greturn, "u", delta*ureturn + ux0, "u0", ux0, "u2", ureturn, "m1", mt0, "m2", mt,"cg", cg, "cgn", cgn, "cgo", cgo
        return - delta*ureturn - ux0

def U2optimized((gx0,gy0),*args):
    """
    U2optimized calculates the sum of the two region two period utilities of regions
    x and y, the other region is myopic. U2optimized passes the relevant parameters
    to U2(gx0,*argx1) and U2(gy0,*argy1) nd returns the sum of these two. See U2 
    for details on the parameters.
    """
    argx,argy = args
    argx[0] = gy0
    argy[0] = gx0
    argx1 = tuple(argx)
    argy1 = tuple(argy)
#    print argx1
#    print argy1
    U2xy = U2(gx0,*argx1) + U2(gy0,*argy1)
    return U2xy


def U3(gx0,*args):
    """
    U3 calculates the two period utility of region x for a given value of
    the other two region's green deployments. Note that U3 depends on the sum of the other region's  indiviual green energy deployments only. 
    U3 is used to calculate the best response of x (i.e.  deployment with max utility) for a given y.\n\n
    Function parameters:\n\n
    t = current period (t+1 = next period)\n
    x = region with a Nash strategy whose best response is to be determined\n
    a,b = the other region with a Nash strategy\n
    c = third region with a myopic strategy\n\n
    gx0 = new green energy of region x in period t (to be determined)\n
    gy0 = sum of new green energy of in regions a and b in period t (given)\n
    Bx = green preference of region x in period t\n
    Bx1 = green preference of region x in period t+1\n
    cnx = levelized cost of new fossil in region x in period t\n
    cnx1 = levelized cost of new fossil in region x in period t+1\n
    cox = levelized cost of old fossil in region x in period t\n
    cox1 = levelized cost of old fossil in region x in period t\n
    Gx0 = total green energy in region x in period t\n
    Gy0 = sum of total green energy in regiona a and b in period t\n
    G0 = total global green energy initially\n
    ex0 = new energy in region x in period t\n
    Ex0 = total energy in region x in period t\n
    Ex1 = total energy in region x in period t+1\n
    M: min( M*Ex0,Ex0 - Gx0) is an upper bound on gx0
    M1: min( M1*Ex1,Ex1 - Gx0 - gx0) is an upper bound on new green energy in region x in period t+1\n
    delta = discount factor of period t+1 wrt period t\n
    detail = boolean flag controlling function return parameters
    """
#    (r2.g[t]+r3.g[t], com.cg0, com.greenpar, r1.B[t], r1.B[t+1], r1.cn[t], r1.co[t], r1.cn[t+1], r1.co[t+1], r1.G[t],\
#    r2.G[t]+r3.G[t], w.G[0], r1.e[t], r1.E[t+1], r1.e[t+1], r1.E[t+2], r1.M[t], r1.M[t+1], (1 - com.discount)**com.years, False)
    gy0,cg0,gpar,Bx,Bx1,cnx,cox,cnx1,cox1,Gx0,Gy0,G0,ex0,Ex0,ex,Ex1,M,M1,delta,detail = args
    mt0 = max(min(M*Ex0,(Ex0 - Gx0)),ex0)
    ux0, cx0 =  0.0, 0.0
    eps = 1e-8
    Gt = Gx0 + Gy0
###    cg = cgx*(Gt/G0)**alpha
    cg = cg2(cg0,G0,Gt,gpar)
#    g1  = Bx0/(cg - cnx) - 1.0
#    g2  = Bx0/(cg - cox) - 1.0
#    cg = cgx*(Gt/G0)^alpha
    cgn = cg - cnx
    cgo = cg - cox
    if cgn > 0.0:
        if mt0 < ex0-eps:
            print mt0,ex0,Gt,Ex0,Gy0,Gx0
            raise ValueError('Max deployment should be greater than next period net demand. Check your assumptions!\n')
    #    if Bx0 <=cg0*(Gt/G0)^alpha - cnx0
        if gx0 <= ex0:
            cx0 = gx0*cgn
            ux0 = Bx*np.log(1.0 + gx0) - cx0 
        elif gx0 > ex0:
            cx0 = (gx0-ex0)*cgo + ex0*cgn
            ux0 = Bx*np.log(1.0 + gx0) - cx0
        else:
            print  gx0, ex0, mt0
            raise ValueError("function U3 1st period: part 1 cgn > 0.0, something wrong! you should never reach here!\n")
    elif cgn <= 0.0 and cgo > 0.0:
        if gx0 <= ex0:
            cx0 = 0.0
            ux0 = Bx*np.log(1.0 + gx0)
        elif gx0 > ex0:
            cx0 = (gx0-ex0)*cgo
            ux0 = Bx*np.log(1.0 + gx0) - cx0
        else:
            print  gx0, ex0, mt0
            raise ValueError("function U3 1st period: part 1 cgo > 0.0, something wrong! you should never reach here!\n")
    elif cgo <= 0.0:
#### This hack prevents the "consumption smoothening"" type behavior that we see in some cases. It forces 
#### ux0 = Bx*np.log(1.0 + gx0) iff gx0 ~= mt0.
        if gx0 >= mt0 -0.5:
            cx0 = 0.0
            ux0 = Bx*np.log(1.0 + gx0)
#            print "SPECIAL CASE", gx0,mt0,ux0
        else:
            cx0 = 0.0
            ux0 = -1e8
#        cx0 = 0.0
#        ux0 = Bx*np.log(1.0 + gx0)
#############################################
    else:
        print  gx0, cg,cox0, cnx0, ex0, mt0
        raise ValueError("function U3 1st period: part 2, something wrong! you should never reach here!\n")
    
#    print "g", gx0, "u", ux0, "e", ex0, "m", mt
    Gx1 = Gx0 + gx0
    Gy1 = Gy0 + gy0
    mt = max(min(M1*Ex1,(Ex1 - Gx1)),ex)
    Gt = Gx1 + Gy1
###    cg = cgx*(Gt/G0)**alpha
    cg = cg2(cg0,G0,Gt,gpar)
    cgn = cg - cnx1
    cgo = cg - cox1
    Bx = Bx1
    g1  = Bx/cgn - 1.0
    g2  = Bx/cgo - 1.0
    
    if mt < ex-0.00001:
        print mt,ex,Gt,Ex,Gy,Gx
        raise ValueError('Max deployment should be greater than next period net demand. Check your assumptions!\n')
    if Bx >= cgn:
        if cgn > 0.0:
            if g1 < ex:
                greturn = g1
                creturn = greturn*cgn 
                ureturn = Bx*np.log(1.0 + greturn) - creturn
            elif g1 > ex and g2 < ex:
                greturn = ex 
                creturn = greturn*cgn
                ureturn = Bx*np.log(1.0 + greturn) - creturn
            elif g2 > ex and g2 <= mt:
                greturn = g2
                creturn = (greturn-ex)*cgo + ex*cgn
                ureturn = Bx*np.log(1.0 + greturn) - creturn
            elif g2 > mt:
                greturn = mt
                creturn = (greturn-ex)*cgo + ex*cgn
                ureturn = Bx*np.log(1.0 + greturn) - creturn
            else:
                raise ValueError("U3 part 2 something wrong! cgn > 0.0: you should never reach here!\n")
#                raise ValueError("nextgG something wrong! cgn > 0.0: you should never reach here!\n")
#            print "cg", cg, "g", greturn, "u", ureturn, "c", creturn, "g1",g1, "g2", g2,"e", ex, "m", mt
#            return greturn, ureturn, creturn, mt
        elif cgn <= 0.0 and cgo > 0.0:
            if  g2 < ex:
                greturn = ex
                creturn = 0.0
                ureturn = Bx*np.log(1.0 + greturn)
            elif g2 > ex and g2 <= mt:
                greturn = g2
                creturn = (greturn-ex)*cgo
                ureturn = Bx*np.log(1.0 + greturn) - creturn
            elif g2 > mt:
                greturn = mt
                creturn = (greturn-ex)*cgo
                ureturn = Bx*np.log(1.0 + greturn) - creturn
            else:
                raise ValueError("U3 part 2 something wrong! cgo > 0.0: you should never reach here!\n")
#                raise ValueError("nextgG something wrong! cgo > 0.0: you should never reach here!\n")
#            print "cg", cg, "g", greturn, "u", ureturn, "c", creturn, "g1",g1, "g2", g2,"e", ex, "m", mt
#            return greturn, ureturn, creturn, mt
        elif cgo <= 0.0:
            greturn = mt
            creturn = 0.0
            ureturn = Bx*np.log(1.0 + greturn)
#            print "cg", cg, "g", greturn, "u", ureturn, "c", creturn, "g1",g1, "g2", g2,"e", ex, "m", mt
#            return greturn, ureturn, creturn, mt
        else:
            raise ValueError("U3 part 2 something wrong! Bx > cgn: you should never reach here!\n")
    elif Bx < cgn:
        greturn = 0.0
        ureturn = 0.0
        creturn = 0.0
#        print "cg", cg, "g", greturn, "u", ureturn, "c", creturn, "g1", g1, "g2", g2,"e", ex, "m", mt
#        return greturn, ureturn, creturn, mt
    else:
        raise ValueError("U2 part 2 something wrong! you should never reach here!\n")
    if detail:
#        print "g", greturn, "u", ureturn, "g1",g1, "g2", g2,"e", ex1, "m", mt
        return  delta*ureturn + ux0, gx0, ux0, cx0, greturn, ureturn, creturn, mt0
    else :
#        print "g", greturn, "u", ureturn, "g1",g1, "g2", g2,"e", ex1, "m", mt,"g0", gx0, "u0", ux0
#        print "gx0", gx0, "g2", greturn, "u", delta*ureturn + ux0, "u0", ux0, "u2", ureturn, "m1", mt0, "m2", mt,"cg", cg, "cgn", cgn, "cgo", cgo
        return - delta*ureturn - ux0

def U3optimized((gx0,gy0,gz0),*args):
    """
    U3optimized calculates the sum of the three region two period utilities of regions
    x, y and z. U3optimized passes the relevant parameters to U3(gx0,*argx), U3(gy0,*argy1) and 
    U3(gz0,*argx) and returns the sum of these two. See U3 for details on the parameters.
    """
    argx,argy,argz = args
    argx[0] = gy0 + gz0
    argy[0] = gx0 + gz0
    argz[0] = gx0 + gy0
    argx1 = tuple(argx)
    argy1 = tuple(argy)
    argz1 = tuple(argz)
#    print argx1
#    print argy1
#    print argz1

    U3_xyz = U3(gx0,*argx1) + U3(gy0,*argy1) + U3(gz0,*argz1)
    return U3_xyz



